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ABSTRACT 

Motivation: Currently there are no curative anticancer drugs, and drug 
resistance is often acquired after drug treatment. One of the reasons is 
that cancers are complex diseases, regulated by multiple signaling 
pathways and cross talks among the pathways. It is expected that 
drug combinations can reduce drug resistance and improve patients' 
outcomes. In clinical practice, the ideal and feasible drug combin- 
ations are combinations of existing Food and Drug Administration- 
approved drugs or bioactive compounds that are already used on 
patients or have entered clinical trials and passed safety tests. 
These drug combinations could directly be used on patients with 
less concern of toxic effects. However, there is so far no effective 
computational approach to search effective drug combinations from 
the enormous number of possibilities. 

Results: In this study, we propose a novel systematic computational 
tool DrugComboRanker to prioritize synergistic drug combinations 
and uncover their mechanisms of action. We first build a drug func- 
tional network based on their genomic profiles, and partition the net- 
work into numerous drug network communities by using a Bayesian 
non-negative matrix factorization approach. As drugs within overlap- 
ping community share common mechanisms of action, we next un- 
cover potential targets of drugs by applying a recommendation 
system on drug communities. We meanwhile build disease-specific 
signaling networks based on patients' genomic profiles and interac- 
tome data. We then identify drug combinations by searching drugs 
whose targets are enriched in the complementary signaling modules 
of the disease signaling network. The novel method was evaluated on 
lung adenocarcinoma and endocrine receptor positive breast cancer, 
and compared with other drug combination approaches. These case 
studies discovered a set of effective drug combinations top ranked in 
our prediction list, and mapped the drug targets on the disease sig- 
naling network to highlight the mechanisms of action of the drug 
combinations. 

Availability and implementation: The program is available on 
request. 

Contact: stwong@tmhs.org 
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1 INTRODUCTION 

Cancers are complex diseases regulated by interactions of mul- 
tiple signaling pathways interacting with each other. Though 
there are some anticancer drugs, they are often not curative 
and associated with acquired drug resistance. It is believed that 
a single drug can target only one disease signaling module, 
whereas alternative signaling pathways will be activated to main- 
tain tumor development. For example, in melanoma, the select- 
ive BRAF inhibitor, PLX4032, can induce tumor responses in 
the majority of patients with BRAF V600E mutant melanoma, 
but the duration of this response is limited owing to the devel- 
opment of acquired resistance, which may be due to alternative 
pathways parallel to the 'BRAF-MEK-ERK' signaling pathway 
(Johannsen et dl., 2013; Nazarian et ai, 2010). To reduce the 
drug resistance effect, and improve cancer treatment outcomes, 
drug combinations are considered as an optimal option. 

The ideal drug combinations expected by clinicians are com- 
binations of Food and Drug Administration (FDA)-approved 
drugs or existing bioactive compounds that have entered clinical 
trials and passed safety tests. These drug combinations could be 
used by patients without toxic side effects. Drug combination 
prediction has been a challenging task in computational biology, 
despite a few studies in the area. The genetic algorithm was, for 
example, proposed to find the optimal combinations of a small 
number of drugs based on the experimental results (Zinner et ai, 
2009). The method, however, is not able to identify optimal drug 
combinations from a large set of drugs, as it requires 
experimental validations in the prediction process. The machine 
learning-based approach was proposed to classify drugs into 
combinations based on molecular and pharmacological features 
of drugs (Zhao et aL, 2011). However, the features of drugs are 
too general to be informative, and the method does not integrate 
the rich genomic profiles of drugs and disease into consideration. 
Combinatorial Drug Assembler (CDA) was developed to iden- 
tify drug combinations based on enrichment analysis of 
genomic profile data. However, the enrichment analysis imple- 
mented in the method is based on differentially expressed genes, 
disease driver genes are not clearly defined, so that drugs target- 
ing the causal dysfunctional signaling modules of disease 
maybe missed out (Lee et aL, 2012). To overcome the afore- 
mentioned problems, we propose a novel computational tool, 
DrugComboRanker, for predicting drug combinations targeting 
multiple signaling modules of cancer-specific networks through 
integrative analysis of genomic profiles of both drugs and 
cancers. 
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2 MATERIALS AND METHODS 

2.1 Overview of DrugComboRanker 

This study seeks to answer a clinical question that which drugs can be 
combined with existing therapy to reduce drug resistance and im- 
prove patients' outcomes. We have thus developed a novel tool, 
DrugComboRanker, to prioritize drug combinations and uncover the 
underlying mechanisms of action. The approach is designed based on 
the following two observations. First, dysfunctional signaling networks 
of diseases are complex and stable, and effective drug combinations can 
inhibit major modules of the networks simultaneously. Second, drugs 
often have multiple active target genes or proteins. Grouping targets of 
drugs with similar transcriptional responses profiles allows better recov- 
ery of their targeted signaling modules in disease than based on indi- 
vidual targets alone. Figure 1 is an overview of the proposed 
DrugComboRanker. The drug functional network is firstly reconstructed 
based on genomic profiling data of drugs that are available in the 
Connectivity Map (CMAP) database. The Bayesian non-negative 
matrix factorization with the /3-divergence (BNMF/3D; Tan and 
Fevotte, 2013) approach is conducted to partition the drug functional 
network into connected functional communities, whose targets indicate 
the inhibited signaling modules of diseases. On the other hand, the dys- 
functional signaling modules of diseases are reconstructed by integrating 
disease genomics data and protein interactome data. Drug combinations 
are then ranked from different drug functional communities, which in- 
hibit disease signaling network modules. 

2.2 Drug functional network reconstruction 

The aim of drug function network reconstruction is to identify drug 
communities that share common transcriptional responses to drug treat- 
ment. We reconstruct, the drug functional network based on the data of 
the CMAP (CMap, build 02; Lamb et al, 2006), which consists of 6100 
gene expression profiles of four cancer cell lines (MCF7, PC3, HL60 and 
SKMEL5) treated by 1309 drugs at different doses. We first combine the 
similarity metric proposed in (Iorio et al, 2010) with the metric defined in 
the STITCH database (Kuhn et al, 2012) as a new drug similarity metric 
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Fig. 1. Overview of the DrugComboRanker for prioritizing disease- 
specific drug combinations 



as follows. First, for each individual drug at each dose, genes are ranked 
based on their fold changes (drug treatment versus control). Then, gene 
rank lists at different doses are merged into one gene rank list by using a 
hierarchical majority voting scheme (Iorio et al, 2010). The merging 
procedure first compares gene rank lists using the Spearman's Footrule 
similarity metric (Diaconis and Graham, 1977), and then the two most 
similar lists are merged repeatedly by using the Borda Merging strategy 
(Lin, 2010). Consequently, gene signatures for individual drugs are cre- 
ated by optimally selecting the top- and bottom-ranked 250 genes, re- 
spectively. This size of these signatures is heuristically determined as 
described in (Iorio et al., 2010). Finally, the gene set enrichment analysis 
(Subramanian et al., 2005) score, So(iJ), is used as the dissimilarity 
metric between drug i and drug j. A threshold T = 0.7848 (the third 
quantile of the empirical probability distribution of the drug dissimilarity 
metric) is given to remove the non-significant dissimilarities as follows: 

fS e (i,A ifS G (UJ)<T 

C, else 

Then we convert the dissimilarity into the normalized similarity score, S> G 
as Sr G (i,j)=l-NORM(Sf G (i,j)) and NORM(Sl G (iJ)) = (Sf G (i,j) 
-min(^(z, .)))/(max(^(z, .)) - wm(Sl G (i, .) where Sf G (i, .) = {Sf G (i, 1), St G 
(i, 2), . . . , S> G (i, n)} denotes all the connection scores linking to the i-th 
node. Here, we empirically set C=10. In addition to the similarity of 
transcriptional responses, a different similarity metric, S s (i,j), was defined 
in STITCH database (Kuhn et al, 2012) by integrating a combined score 
of structure similarity, experimental similarity and text mining similarity 
score. The structure similarity is defined by the Tanimoto 2D chemical 
similarity scores (Martin et al, 2002). The experimental similarity is cal- 
culated by the Pearson correlation of the activity patterns of the com- 
pounds based on NCI60 screens. The text mining similarity is computed 
by mining the curated database, such as OMIM and MEDLINE, using a 
co-occurrence scheme and a natural language processing approach 
(Jensen et al, 2006; Saric et al, 2006). In this study, we combine the 
two similarity metrics as S(i,j) = S s (i,j) + Sf' G (i,j). The rational to use 
the combined similarity is to generate better drug community structure 
when applying the clustering algorithm on the drug functional network. 
Finally, a drug functional network is reconstructed with 1308 nodes and 
51 164 edges. Figure 2 shows an example of the reconstructed drug func- 
tional network community. 

2.3 Drug community discovery by using Bayesian non- 
negative matrix factorization with ^-divergence 

Next, we partition the reconstructed drug functional network into com- 
munities (modules) within which drugs share common mechanisms of 
action (targeting the same or related signaling modules), as illustrated 
in Figure 2. We use the BNMF/3D (Tan and Fevotte, 2013) for the net- 
work partition. The algorithm has shown good partition results on vari- 
ous classes of real-world signals such as audio and images (Tan and 
Fevotte, 2013), Specifically, the algorithm can (i) enable soft partitioning 
(one drug can be partitioned into multiple communities with different 
membership scores); (ii) reduce noise influence (suppresses noise in the 
high dimensional space); and (iii) allow multi-scale partition. To the best 
of our knowledge, the algorithm is for the first time used in drug func- 
tional network partition. 

Specifically, given a weighted undirected network with N nodes, let 
V e [R+ xN denote its adjacent matrix. BNMF/3D facto rizes the matrix V 
into two low-rank matrices: 

\^\ = WH (1) 

where W and H are non-negative matrices with dimensions N x K and 
K x N, and the elements of V, V, W and H are denoted by vy, hj,^ik and 
hkj, respectively. The /c-th column of Wis denoted by Wk e [R+ xK and the 
k-th row of //is denoted by ht 6 [R+ xN . The Wis normalized to W^, that 
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Fig. 2. Snapshots of the drug functional network (top) and an example of 
the reconstructed drug network community (bottom) 



is, each row of Wn is subject to 



7=1 



where wriy is the element of Wn and it quantifies the membership of node 
i with respect to community j. To infer the appropriate model order K, 
we use a Bayesian paradigm for non-negative matrix factorization 
by placing automatic relevance determination priors with scale hyper- 
parameters A = {kk}£=i on the variables w// c and hq. In this model, the 
distribution of X k is parameterized by fixed parameters a and b, and the 
fixed parameter ft decides the distance measure between the observed 
interactions V and the expected interactions V. Under these assumptions, 
the posterior density function can be obtained as 

_p(V\W,H)p(W\X)p(H\l)p(X) 



p(W,H,k \ V)~- 



P(V) 



Maximizing the posterior density is equivalent to minimizing the negative 
log posterior, which can be regarded as a loss function Cmap(W, H, X) as 



C MAP (W, H,X)= - log p(W, H,k\V) 

= - \ogp(V | W, H) - log p(W\k) - \og(H\k) - log p(k) 

where log p(V\ W, H) is the log-likelihood. 
The generalized ^-divergence is defined by 



(4) 
(5) 



D p (x\y)± 



xy H 



1) p p-1 



,PeM\{0, 1} 



xlog x + y, i 

x , x „ n 
--log- -1,0 

y y 



= o 



The ^-divergence can be regarded as a minus log-likelihood for the 
Tweedie distribution and its probability density function is given by 



fix, n,<i>,P) = h(x, 0)exp | \ (— '- 



xn H 



1 



(7) 



where h(x, <j>) is the base measure function, /x is the mean, 0 is the dis- 
persion parameter and ft is the shape parameter. Assuming that v# is 
generated from the Tweedie distribution, the log-likelihood function 
can be given by 



-log p( V | W, H) = i Dp( V | WH) + C 



(8) 



To insure W and H are non-negative, the Half-Normal priors are as 
signed on them, 

Pim \X k ) = HN(w ik \X k ) 



where 



HAf(x \X) 



p{h kj \X k )=HAf(h kj \X k ) 
2 \i 



and place an inverse Gamma priors on each X^, 



p(X k ; a, b)-- 



r(«) 



'exp(--) 



(9) 
(10) 

(11) 
(12) 



Then, according to Equation (5), the objective function Cmap(W, H, X) 
can be given as 



C MAP i W,H,)=^Dp(V\WH) + ^Y k (^ w l+\ h l + 



k=\ 

+ (N+a+l)\ogX k + C 



(13) 



To minimize C M Api W, H, X) with respect to W, H and k, we adopt the 
strategy in (Tan and Fevotte, 2013) by introducing a local majorization- 
minimization algorithm with efficient multiplicative updates. Finally, we 
give the overlapping community detection algorithm as follows: 

Step 1: Initialize Wk e U^_ xK andhk £ [R+ xN to random non-negative 
values. 



Step 2: Update W, H and X k by 

W T HWH)^- 2) ■ V) 



H=H 



W=W 



W T iWH) <p - l) + 4>/repmat(JL, 1, AO 



((WH) 



■(fi-2) 



V)H l 



(W^-^H 1 +(f)/repmatik, \,N) 
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l/Q-fi,P<2 
l/(fi-l),p>2 



Step 3: Repeat Step 2 until max k =\,2,...,K\{K ew ^° k ld )/^ k d \ < £ 
Step 4: Normalize W to Wn, then the number of non-zero columns 
k* of Wn is the number of clusters. Assign each node to the k* 
clusters according to Wn- In the above algorithm, X- Y denotes 
element-by-element multiplication of X and Y; j denotes element- 
by-element division of X and Y; and X y denotes raising each element 
of X to the y th power. In addition, repmat(k, 1 , AO denotes the Kx N 
matrix with each column being the vector^. Using the BNMF/3D 
approach, we partition the drug functional network into a set of 
connected network modules (Fig. 2), within which drugs share 
common targets or related signaling mechanisms. 



2.4 Drug combination discovery based on target network 
analysis 

The novel drug combination approach consists of the following three 
major components. 

(a) Disease specific signaling network reconstruction 

Several approaches (Barrenas et al, 2012; Chuang et al, 2007; Ideker 
et al, 2002) have been proposed to reconstruct signaling networks of 
diseases based on transcritpome and interactome data. In this study, 
the approach proposed in (Barrenas et al, 2012) is used. The integrated 
protein-protein interactions from BioGRID (Stark et al., 2006) and the 
manually curated human cancer signaling from (Awan et al, 2007; Cui 
et al, 2007; Li et al, 2012; Newman et al, 2013) (available at http://www. 
bri.nrc.ca/wang/) are clustered into functional protein-protein modules, 
and each module is tested for enrichment (Fisher's exact test, P<0.01) of 
the differentially expressed genes of the gene expression profile of a dis- 
ease. The enriched modules are then considered as disease susceptibility 
modules, from which the highly interconnected genes are identified as the 
disease-specific signaling network. Figure 3 shows the reconstructed sig- 
naling network of lung adenocarcinoma. 

(b) Functional drug target prediction using network-based 
recommendation 

Drugs often have multiple targets and affect distinct signaling modules, 
but only parts of them are known for given drugs. The drug communities 
embed targeting signaling modules (functional targets instead of physical 
targets) of drugs. To uncover the targeting signaling modules of drugs, 
we propose a network-based recommendation approach as follows. Let 
D = {d\, da, . . • , d m ] denote drugs in a given community, and T= 
{t\, t2, ■ ■ ■ , t n ) be the known drug targets. The drug-target interaction 
network can be described as a bipartite graph G(D, T, E). The E indicates 
the known drug and target associations. This drug-target network can be 
represented by an adjacent matrix A= {%} nxw > where aji is the weight 
that quantifies the association between d\ and tj. Then the novel network- 
based algorithm is designed based on a bipartite network projection tech- 
nique (Zhou et al, 2010) as follows. 



R = FA 



(14) 



where R= {ry\ is the recommendation score (the functional associ- 
ation possibility between drug i and target j). The F= {/ly} indicates 
the transition matrix from drug i to drug j and is defined as: 



r(i,j)f-;k(xi) 



E 



(15) 




o U O O ' 




Fig. 3. The reconstructed signaling network of lung adenocarcinoma. 
Different node colors indicate different modules, and the node size indi- 
cates the degree of nodes 



where r(i,j) = k(t i ) 1 ~ x k(t J ) and k(x) is the degree of the x node in the 
bipartite network. Targets with recommendation scores greater than 
given threshold, 0.1, are kept as the active functional targets of given 
drugs. 

(c) Disease-specific drug combination discovery 

With the constructed disease signaling network and the predicted 
drug targets, drug combinations are then prioritized by combining the 
following synergistic scores. Given two candidate drugs, d t and dj, 
from different clusters, suppose d t e C& and dj e C/ 7 , and 7& = 
{tki, tki, ■ ■ ■ , tkm} denote the targets of d t in Q, and Tj x = 
{thi, hi, • • • , hn} denote the targets of dj in C/ 7 . The first synergistic 
score is defined as follows. 



J2iCS(t ki )ex.p(-D {t, d , T h }/n 2 ) 



(16) 



where CSfad is the centrality score of target tu in the reconstructed 
disease signaling network, and it is an additive of betweenness (Bn), 
closeness (Cn) (Brandes and Fleischer, 2005) and PageRank (Pr) score 
(Page et al., 1999) of protein t ki , that is 

CS(t ki ) = Bn(t ki ) + Cn(t ki ) + Pr(t ki ) (17) 

These are three different but correlated centrality measurements, and the 
reason of combing them is to get a robust centrality score. The min D 
{t k i, Th] is the minimum shortest path from t ki to T^. The first synergistic 
score, Si(i,j), prefers drug combinations, whose targets are in the center 
(hubs) of disease signaling network and closely connected. 
The second synergistic score is defined as 



S 2 (i,j) = 



^2 Sim(tki, hj) 
(m + n)(m + n — 1) 



(18) 
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where Sim{tki,thj) is the semantic similarity of gene ontology (GO) anno- 
tations of tki and tjy (Couto et al, 2007), which is computed based on 
the overlap of the GO terms that are associated with and and is 
defined as 



Sim(t ku t hj ) 



21og 2 m<2x{X^)} 



(\o% 2 p(GO ki ) + \o% 2 p(GO hj j) 



(19) 



Where GOki is the GO term that associated with t^u and A is a GO term 
that is an ancestor of both GOki and GOhj, and 

Freq(GO ki ) 



p(GO ki ) = 



MaxFreq 



(20) 



Freq{GOki) is the frequency of GO term GOki occurring in GO annota- 
tions, which are taken from GO database. MaxFreq is the maximum 
occurrences frequency of GO terms that are associated with all the targets 
and the predicted drug targets in the GO annotations. The second syn- 
ergistic score, SjiiJ), prefers drug combinations that block genes with 
similar functions, e.g. cell proliferation. 

Finally, the synergistic score of drug d t and dj is given by 

S y (i,j) = (S^iJ) + S 2 (i,j)) ■ d E (iJ) (21) 



where 



d E (i,j) = max{s G (i,j), ^S G (i,j) ■ SJijj, S G (iJ) + S s (i,j)} 



(22) 



where S G (i,j) reflects the distance of their expression pattern. In sum- 
mary, drug combinations targeting on the disease-specific signaling 
network, with similar functions, through alternative targets are 
prioritized. 



3 RESULTS 

We have applied the BNMF/3D algorithm to cluster the drug- 
drug network into overlapping drug communities. Table 1 lists 
the parameters in the analysis. We set a as one of 
{10,50,100,150,200,250,300,350,400,450}, and b to be equal to 
a. When a and b are set to be 450, the BNMF^D algorithm 
will converge more quickly. Under this setting of a and b, we 
test the drug network community reconstruction results with dif- 
ferent K and ft as K = {100,110,120,130,140,150,160,170} and 
/? = {1, 2}. When K'\s 150 and f$ is 1, we get better drug commu- 
nity structures in terms of the following defined community qual- 
ity measurements. We assign the number from 1 to 150 to each 
drug community, and obtain the drug-target networks by con- 
necting the drugs in the same community through their targets. 
The connected drug-target networks provide us a useful way to 
drug discovery. 



Table 1. The values of parameters in BNMF/JD algorithm 



Parameters 


Values 


K 


150 


P 


1 


a = b 


450 


£ 


lO" 5 



3.1 Quality evaluation of predicted drug community 

To evaluate the proposed approach, we compared it with three 
widely used network partition approaches: Markov cluster 
(MCL; Enright et al, 2002), affinity propagation (AP; Frey 
and Dueck, 2007) and CluterONE (Nepusz et al, 2012). To 
evaluate the quality of the predicted drug community, the fol- 
lowing three metrics are defined. The average enrichment score, 

siNi p 

S N , is defined as S N = T^fc 1 — , where s t = jj-, P t is the number of 

drug pairs sharing targets in the z'-th drug community, N t is the 
number of drugs of that drug community and N is the number of 
drug community. Higher S N value indicates more drugs in the 
same community sharing common targets. The metric, Sr, is 
defined as the relative number of well-connected drug commu- 
nities Sr = MA^M^ where |.| is the number of elements, and d is 
a given threshold. The geometric accuracy metric, Sga, is defined 
as the balance of the two metrics: Sga = V^V • Sr. The compari- 
son result shown in Figure 4 shows that the proposed method 
outperforms the other three methods in terms of all three metrics 
with d=0.6. 



3.2 Drug combination prediction evaluation 

To evaluate the prediction capacity of DrugComboRanker, 
we apply it on the lung adenocarcinoma and endocrine 
receptor (ER) positive breast cancer. The 107 lung 
adenocarcinoma patient samples (58 lung tumor tissues and 49 
normal lung tissues, GDS3257) were collected to reconstruct 
the lung adenocarcinoma- specific signaling network, as 
shown in Figure 3. This signaling network consists of 
905 genes and 1 1 282 interactions. Based on the signaling 
network, we predict effective drugs that have synergistic effects 
combined with three drugs (Gefitinib, Paclitaxel and LY-294002) 
given by our clinicians. Level 3 gene expression data 
of 222 ER positive breast invasive carcinoma (BRCA) 
samples and 61 normal breast samples were collected 
from TCGA. The synergistic drug combination with two 
FDA- approved drugs, Tamoxifen and Letrozole, are further 




Fig. 4. The comparison results of the proposed method (referred to as 
BNMF PD) with AP, ClusterONE and MCL according to the quality of 
the predicted drug communities in terms of S N , S R and S GA 
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Table 2. Synergistic alternative drugs combining with Gefitinib 



Drug combination 


Synergistic 


Rank 


Literature 


(Community number) 


score 




evidence 


Gefitinib/Paclitaxel (63/55) 


2.905 


2 


PMID:19596955 








PMID: 14990633 


Gefitinib/Celecoxib (63/48) 


2.804 


3 


PMID:18379355 








PMID: 16914589 


Gefitinib/Genistein (63/102) 


2.765 


4 


PMID:22160570 


Gefitinib/Fulvestrant (63/55) 


2.529 


9 


PMID:24268810 


Gefitinib/Irinotecan (63/66) 


2.468 


11 


PMID:21915126 








PMID:16713012 


Gefitinib/Vorinostat (63/102) 


2.464 


12 


PMID:21271222 


Gefitinib/Lovastatin (63/34) 


2.102 


27 


PMID:19760159 


Gefitinib/Rosiglitazome (63/13) 


2.023 


32 


PMID:168386327 


Gefitinib/MS-275(63/102) 


2.007 


34 


PMID:16424029 


Table 3. Synergistic alternative drugs combining with LY-294002 


Drug combination 


Synergistic 


Rank 


Literature 


(Community number) 


score 




evidence 


LY-294002/Quercetin (55/140) 


2.987 


2 


PMID: 14688022 


LY-294002/SB-202190(55/5) 


2.936 


3 


PMID: 161 15952 


LY-294002/Rapamycin(55/54) 


2.448 


33 


PMID:22355375 


Table 4. Synergistic alternative drugs combining with Paclitaxel 


Drug combination (Community 


Synergistic 


Rank 


Literature 


number) 


score 




evidence 


Paclitaxel/Irinotecan(5 5/66) 


3.096 


2 


PMID: 1625 1879 


Paclitaxel/Gefitinib(55/63) 


2.905 


8 


PMID: 19596955 








PMID: 14990633 


Paclitaxel/Vorinostat(55/102) 


2.901 


9 


PMID: 1962 1389 


Paclitaxel/Rapamycin(55/54) 


2.888 


11 


PMID:22896668 


Paclitaxel/Vinblastine(55/1 37) 


2.681 


17 


PMID:16193638 


Paclitaxel/Ifosfamide(55/90) 


2.377 


25 


PMID:8740798 


Paclitaxel/SuraminSodium (55/3) 


2.189 


37 


PMID: 11 507065 


Paclitaxel/Resveratrol(55/48) 


2.060 


49 


PMID: 146667 16 



predicted for ER-positive BRCA based on the reconstructed ER- 
positive breast cancer signaling network. 

To validate the prediction results, we searched the literature 
evidence of the top 50 ranked combinations. Tables 2-4 show the 
literature evidence of the effective drug combinations in our top 
50 lists for lung adenocarcinoma. Surprisingly, 19 different drug 
combinations have been reported to be synergistic combining 
with Gefitinib, Paclitaxel and LY-294002 in non-small-cell 
lung cancer. Tables 5 and 6 show the literature evidence of ef- 
fective drug combinations in the top 50 lists for ER-positive 
breast cancer. Also 14 different drug combinations have been 



Table 5. Synergistic alternative drugs combining with Tamoxifen 



Drug combination 


Synergistic 


Rank 


Literature 


(Community number) 


score 




evidence 


Tamoxif en /Tretinoin (49 / 5 5) 


2.913 


8 


PMID:9533531 


Tamoxifen/Genistein(49 / 1 02) 


2.910 


9 


PMID: 17295235 


Tamoxif en/Quercetin(49/140) 


2.831 


12 


PMID:20804812 


Tamoxif en/Verapamil(49/66) 


2.776 


19 


PMID:8 185686 


Tamoxif en/Mifepristone(49/66) 


2.721 


21 


PMID:9879777 


Tamoxifen/Celecoxib(49/48) 


2.437 


29 


PMID:23731702 


Tamoxif en/Resveratro 1(49/48) 


2.416 


31 


PMID:23896596 


Tamoxif en/Letrozole(49/92) 


2.299 


36 


PMID: 15026471 


Tamoxif en/Exemestane(49/51) 


2.136 


42 


PMID:2 1247627 


Tamoxifen/LY-294002(49/55) 


2.035 


46 


PMID:12479367 



Table 6. Synergistic alternative drugs combining with Letrozole 



Drug combination 


Synergistic 


Rank 


Literature 


(Community number) 


score 




evidence 


Letrozole/Fulvestrant(92/55) 


2.457 


8 


PMID: 15958593 


Letrozole/Tamoxifen(92/49) 


2.300 


10 


PMID: 15026471 


Letrozole/Metf ormin(92 /56) 


2.148 


15 


NCT01 589367 


Letrozole/Celecoxib(92/48) 


2.128 


16 


PMID:19254941 


Letrozole/Imatinib(92/55) 


1.765 


40 


NCT00338728 


reported to be synergistic 


combining 


with 


Tamoxifen and 



Letrozole. These results show the strong drug combination pre- 
diction capacity of DrugComboRanker. 

To make the evaluation fair and sound, we further compared 
the predicted results with CDA, and a random combination 
method (RCM), which randomly picks up 50 drugs from the 
available drug lists to combine with the designated drugs. As 
for CDA, we picked the top 50 drug combinations with the 
designated drugs. Figures 5 and 6 show the comparison results 
on the lung adenocarcinoma and ER-positive breast cancer in 
terms of literature supports of those top-ranked 50 drug com- 
binations, respectively (The numbers on the bars of RCM are 
standard deviations). For RCM, we repeated the random selec- 
tion 100 times; for each simulation, we checked the literature 
evidence (In total, we checked all the 5000 random combin- 
ations). As can be seen, the proposed approach outperforms 
the CDA and random selection significantly. 

The predicted drug targets in the disease-specific network 
could indicate the molecular mechanism of synergistic drug com- 
binations. Here, we map the responsive genes of Gefitinib, 
Paclitaxel, Vorinostat, LY-294002 and Quercetin to the lung 
adenocarcinoma- specific signaling network to capture the dis- 
tinct synergistic responses induced by three agent combinations, 
Gefitinib and Paclitaxel, LY-294002 and Quercetin, Gefitinib 
and Vorinostat. As shown in Figure 7, Gefitinib and Paclitaxel 
combinations can affect the EGFR signaling pathway (endothe- 
lial cell proliferation), TP53 signaling pathways, as well as 
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in terms of the number of literature supports of the top-ranked 50 drug 
combinations of lung adenocarcinoma with designated drugs, Gefitinib, 

Tamoxifen-i- Letrozole+ 

Paclitaxel and LY-294002 

Fig. 6. The comparison results of DrugComboRanker, CDA and RCM 
in terms of the number of literature supports of the top-ranked 50 drug 
combinations of ER-positive breast cancer with designated drugs, 
Tamoxifen and Letrozole 




Fig. 7. Drug targets mapped on the disease signaling network. Red and green are the drug targets of Gefitinib and Paclitaxel, respectively 
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Fig. 8. Drug targets mapped on the disease signaling network. Red and green are the drug targets of LY-294002 and Quercetin, respectively; blue nodes 
are the weak effected targets of both drugs 




Fig. 9. Drug targets mapped on the disease signaling network. Red and green are the drug targets of Gefitinib and Vorinostat, respectively 
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biological processes, such as cell cycle, apoptosis and the hub 
genes, i.e. EGFR, TP53, SRC, FOS, JUN. Distinctly the LY- 
294002 and Quercetin combination affects the alternative 
EGFR, PI3K-AKT and JAK-STAT3 pathways, as can be seen 
in Figure 8. In addition, the drug combination, Gefitinib and 
Celecoxib, targets the EGFR and COX-2 signaling pathways, 
respectively. The Gefitinib and Celecoxib have distinct transcrip- 
tional responses that indicate EGFR and COX-2 signaling path- 
ways are complementary, and have cross talks. Another example 
is Gefitinib and Vorinostat combination. Vorinostat is a Histone 
deacetylases inhibitor, as shown in Figure 9, which interacts with 
CTNNB1, and CTNNB1 interacts with E-cadherin, ERBB2 and 
EGFR, whereas Gefitinib targets on EGFR. Thus, this combin- 
ation forms a double inhibition on growth factors. 

4 DISCUSSION AND CONCLUSION 

Drug combinations can provide an optimal cancer therapeutic 
treatment by overcoming acquired drug resistance. Clinicians 
usually choose drug combinations manually based on their ex- 
perience and expertise. However, this expert-picking-based 
method cannot be scaled up for discovering drug combinations 
from a large set of drugs. The availability of genomic profiles of 
drugs and patients are informative for uncovering the regulatory 
signaling networks of cancers and the mechanisms of action of 
drugs. However, the computational tools are the bottleneck to 
interpret and convert the big genomics data into discovery of 
drug combinations. In this study, we propose a novel computa- 
tional tool, DrugComboRanker, to prioritize drug combin- 
ations for specific cancers by selecting combinations targeting 
the alternative and complementary signaling modules of disease. 
We first predict the functional targets of drugs via the network 
community analysis of drug functional networks based on their 
genomic profiles. We then define the synergistic score of drug 
combinations that block key disease genes and their potential 
cross talks. The new method can also provide insights into mech- 
anism of actions of drug combinations by mapping the pre- 
dicted drug targets on the disease signaling network. In our 
case studies on lung adenocarcinoma and ER-positive breast 
cancer cases, a set of top-ranked drug combinations, with similar 
and distinct mechanisms of action, top-ranked in our predic- 
tion list have been reported to be effective in anticancer treat- 
ment. Moreover, a set of novel drug combinations are also 
identified with a better anticancer treatment outcome. In conclu- 
sion, the DrugComboRanker has the potential to be a widely 
used software package for predicting new disease- specific drug 
combinations. 
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